The Thermodynamical Limit of the Lipkin-Meshkov-Glick Model 
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A method based on the analyzis of the Majorana polynomial roots is introduced to compute the 
spectrum of the Lipkin-Meshkov-Glick model in the thermodynamical limit. A rich structure made 
of four qualitatively different regions is revealed in the parameter space whereas the ground state 
study only distinguishes between two phases. 
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The Lipkin-Meshkov-Glick (LMG) model introduced 
in 1965 to describe shape phase transition in nuclei [lj 
has, since then, been proposed to describe many sys- 
tems ranging from interacting spin systems Q to Bose- 
Einstein condensates 0] or magnetic molecules such as 
Mni2 acetate This ubiquity is due to its mapping 
onto a single-particle evolving in a double-well potential 
0,0] or onto an interacting two- level boson system. More 
recently, this model has also been used to investigate the 
relationship between entanglement and quantum phase 
transitions 0, & S, • 

The LMG model is known to be exactly solvable 



matrices: 



id, E3, [13]. However, getting the solution requires to 
solve Bethe-like equations, a task which, in the present 
context, is more costly than exact diagonalization. A 
complete description of the spectrum thus requires to 
develop alternative routes. Though the low-energy spec- 
trum has been studied in details via different meth- 
ods (variational 0] , bosonization 0, E3] , coherent states 
[TBI), the richness of the full spectrum has been investi- 
gated only lately by means of numerical diagonalizations 
flil E3] • These latter studies suggest the existence of sin- 
gular points in the density of states as well as a nontrivial 
level spacing distribution. 

In this paper, we shed light on these issues by exactly 
computing the spectrum of the LMG model. The pro- 
posed method relies on the determination of the Majo- 
rana polynomial roots associated to the eigenstates of the 
Hamiltonian. This polynomial is built within a coherent 
state formalism which is well-suited to such a system. 
Within this framework, the spectrum is encoded in a lin- 
ear differential equation which is solved in the thermo- 
dynamical limit. This allows us to exactly compute the 
density of states in the whole parameter range and to 
locate its singularities. Four distinct regions arise with 
qualitatively different properties. In particular, we find a 
parameter regime for which the density of states has no 
thermodynamical limit. 

The LMG model describes a set of N spins | mutually 
interacting through a XK-like Hamiltonian and coupled 
to an external transverse magnetic field h. This Hamil- 
tonian H can thus be expressed in terms of the total spin 
operators S a = YliLi where the ct q 's are the Pauli 



H = --( lx Sl 



(1) 



In the following, we only consider the maximum spin 
sector s — N/2 with N even. Given the symmetry of 
the spectrum of H, we focus on the parameter range 
h ^ 0; s$ |7„| < j x . Note also that [H, S 2 ] = and 



[H, 



,i7r(S z — s) 



] = (spin- flip symmetry). Denoting by 



{\s,m}} the standard eigenbasis of {S 2 ,^}, this latter 
symmetry implies that odd and even states decouple. In 
the thermodynamical limit, both subspaces are isospec- 
tral so that we limit the following analyzis to the sector 
m even for which one has exactly (s + 1) eigenstates. 

In the spin coherent states basis |l8j|, with non- 
normalized states | a) = e a5+ |s,— s), any state |$) is 
represented by its Majorana polynomial [19( defined as 

(s,m\y)a m+s , 
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where d $J 2s is its degree. 

The standard representation of the spin operators 
(S± = S x ± iS y ) in the coherent states basis 

S+ = 2sa~a 2 d a , (3) 
S- = B a , (4) 
S z = -s + ad a , (5) 

allows to map the Schrodinger equation H\^)= E\^f) 
onto the following linear differential equation 0, H3] 
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Except for trivial values of the parameters, the degree of 
VP for an eigenstate of H in the sector we considered is 
d = 2s. At this step, the spectrum could be analyzed by 
mapping Eq. ([5]) onto a Schrodinger equation describing 
a particle in an effective one-dimensional potential 0, Q . 
Then, a semi-classical treatment would, in principle, al- 
low one to obtain the density of states in the thermo- 
dynamical limit, as shown in Ref. 0] for the low-energy 
spectrum in the region ^ y = 0, j x > 0. Unfortunately, for 
arbitrary values of the parameters, the effective potential 
becomes tricky Q and such an approach is therefore dif- 
ficult to follow. Here, we propose an alternative route by 
first converting the linear second-order differential equa- 
tion for 'J ([6]) into a first-order differential equation for 
its loga rithmic derivative. More precisely, the function 
G defined as 



2s 

^) = |^log*(a) = i£^ (10) 

k—1 



satisfies the following Riccati-like equation 
'G'(a) 



2s 



+ G>) 



+ P 1 (a)G(a) + P {a)=e. (11) 



The density of states is then given by analyzing the poles 
of G, i.e., the roots of the Majorana polynomial \F In- 
deed, the cornerstone of this study is that, for this model, 
the afc's are spread over two curves Cq and C\ in the com- 
plex plane which depend on the energy. In addition, the 
n-th excited state of H has 2n poles on C\ and 2(s— n) on 
Co (thus defining both curves). This remarkable property 
stems from the oscillation theorem which indexes the ex- 
cited states for a particle in the effective one-dimensional 
potential (discussed above) by the number of wavefunc- 
tion nodes. To illustrate this repartition of the poles 
which is likely related to the integrability of the model 
[2ll ] . we display in Fig. [1] several typical states in the 
Majorana sphere representation [l9| . This representation 
generalizes the celebrated Bloch sphere used for spin I 
states and proceeds as follows. For a given polynomial 
with d roots, we first complement it with (2s — d) roots 
at infinity in the complex plane. Next, the resulting set 
of 2s points is sent onto the unit sphere by an inverse 
stereographic map. For instance, within this mapping, 
the basis state |s, m) is represented by (s — m) points on 
one pole and (s + m) points on the opposite pole. 

The location of the poles of G explained above pro- 
vides a straightforward relation between the normalized 
integrated density of states N G [0, 1] and the number of 
poles lying in C\ . One indeed simply has 
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where C\ is a contour that surrounds C\ and oriented such 
that J\f > 0. In this equation, G is built from the n-th 




FIG. 1: Upper part: representation of the poles of G on 
the Majorana sphere (blue dots) for three typical eigenstates 
computed for h = 1, 7% = 5, = —3 and s = 20 (zone III in 
Fig. [2J. Black lines correspond to the Go branch cuts Co and 
C x . 

Lower part: Numerical (black staircase curve s = 20) ver- 
sus analytical (red line s — 00 ) integrated density of states. 
Black dots indicates the singularity of the density of states 
N& n (-h) and N^ n (h) (Eqs. JT5J and |gT} respectively) in 
the thermodynamical limit. 



excited state and the dependence of AT with e is given 
via Eq. (fTTj) . Unfortunately, one cannot solve Eq. (fTT|) 
exactly at finite s, which would give a complete explicit 
solution to our problem. However, one can easily solve it 
perturbatively in 1/s. 

Therefore, let us assume that G, and e, can be ex- 
panded in the following form 

G, 



G 



E 

ten 



(13) 



At leading order s°, Eq. (fTTj) becomes a second-order 
polynomial equation for Go whose solutions simply read 



G ± (a) 



where 



t [a 2 ( lv - lx ) + 7x + 7y + 2/i] ± ^/2Q{a) 
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(14) 



Q(a) = (jy--fx){h + eo)a 4 + 2[h 2 +j x -f y + e (jx+7y)]a: 2 
+ (7x-7»)(fc-eo). (15) 

The four roots of Q are the branch points of G which de- 
fine the limit of the curves Co and C\ . A close analyzis of 
these branch cuts, in the parameter space, then leads to 
the integrated density of states in the thermodynamical 
limit which reads 

Urn JV(e) = Afo(eo) = ^ / [G+(a) - Gq (a)] da. 

(16) 
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This quantity can be expressed in terms of lengthy ex- 
pressions involving elliptic integrals. We obtained these 
expressions in the whole parameter space and will give 
them explicitely in a forthcoming publication [22| . In 
the following, we only discuss qualitatively the various 
regions that must be distinguished, and analyze them by 
means of the density of states po(eo) = <9 £o Ao(£o)- 

But, first of all, let us remind the reader that if one 
only considers the properties of the ground state, which 
define the zero-tempcrature phase diagram, only two 
phases must be distinguished. For h > j x (symmetric 
phase), the ground state is unique and lim^oo (S z ) /s = 1 
whereas for h < j x (broken phase), the ground state 
is two-fold degenerate and lim^^ (S z )/s — h/j x . The 
quantum phase transition at h = ^ x is second-order and 
characterized by mean-field critical exponents and 
nontrivial finite-size scaling behavior [ij, [23[ . An impor- 
tant result of our study is that, when considering the 
full spectrum, four different zones arise instead of two. 
These regions, described below, are characterized by dif- 
ferent singular behaviors of the density of states (see Fig. 
(2) as already noticed in a numerical study of the special 
case 7z = -j y 

• Zone I: \~f y \ < y x < h. In this sector, the density of 
states po is a smooth function of — h ^ Sq ^ h as can be 
seen in Fig. [2] The distribution of Majorana polynomial 
roots for the eigenstates is similar to that displayed in 
Fig. QJ6). In the complex plane, Co and C\ lie in the 
imaginary and real axes respectively. 

• Zone II: \ r y y \ < h < j x . In this region, two distinct 
branches must be distinguished: 



11(a): 



27. 



$J Eq —h. Cq coincides with the 



whole imaginary axis while C\ is made of two discon- 
nected segments in the real axis as depicted in FigOJa); 
— 11(6): —h ^ £o ^ h. Co and C\ are the same as in I. 
These two branches of the density of states diverge 
at £o = —h where the elliptic integrals involved in the 
expression of Ao can be rccasted in the simple following 
form 
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+ \fWy + ^{lx-h){h- ly ), (19) 
±(V^h - VWy) + \Jilx- h)(h - 7 „). (20) 



• Zone III: h < —~f y < j x . In this zone, there are three 
different branches: 

— 111(a): ^ £ o ^ — h. Cq and C\ are the same 

as in 11(a); 
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FIG. 2: Phase diagram in the j x , y y plane at fixed h > and 
typical density of states for ( r fx, r ^y, h) equal to I: (1/2, 1/3, 1), 
IL (2,1/2,1), III: (5,-3,1), IV: (5,3,1). 



— 111(6): —h ^ £o h. Co and C\ are the same as in I; 

h 2 +*y 2 

— III(c): h < e Q < ^ JL - Co is made of two discon- 
nected segments in the imaginary axis while C\ coincides 
with the whole real axis as depicted on the Majorana 
sphere in FigfTJc). 

In this zone, the density of states has two singularities 
at Eo = ±/i. Their position in the spectrum is given by 
No U (-h) =Af ll (-h) [see Eq. JTHJ] and 



A+ tan" 



A + 

-l A h 



A h tan 



(21) 

For 7 X = — 7 y , the density of states is symmetric with 
respect to £o = 0, and the above expression gives the 
exact location, in the thermodynamical limit, of the so- 
called exceptional point observed in Ref. [161 ] . 

• Zone IV: h < j y < j x - This part of the phase dia- 
gram is the most complex one. The density of states is, 
as in zone III, made of three different branches: 
h 2 +1 l ^ _ ^ h 2 +j 2 



IV(d): 

IV(e): 
IV(6): 



27* 
h 2 +7? 



27y 



27y 

-h ^ £o ^ h. 

However, the structure of the C\ curve, in each case, 
is complex but Co remains simple so that the integral in 
Eq. p6p can still be computed and reveals two special 

/i 2 +7 2 

points. The first one occurs at Eq = — - 
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In conclusion, we would like to emphasize that the 
present approach can be extended to other similar models 
with higher-order interaction terms where the mapping 
onto a particle in a one-dimensional potential fails. Fur- 
ther, one can go beyond the thermodynamical limit and 
extract the finite-s corrections which could be crucial for 
some observables [22j]. 

We are grateful to C. Aslangul, S. Dusuel, N. Gromov, 
J.-M. Maillard and P. Vieira for fruitful and stimulating 
discussions. PR was partially supported by FCT and EU 
FEDER through POCTI and POCI, namely via Quant- 
Log POCI/MAT/55796/2004 Project of CLC-DM-IST, 
SQIG-IT and grant SFRH/BD/16182/2004/2ZB5. 



FIG. 3: Gap between two consecutive levels as a function of 
the energy in region IV for j x — 15, y y = 10 and h = 1. In 
the central region, the red line is the average gap Ao in the 
thermodynamical limit. The branches below and above this 
line are respectively A' - ^ and A^ (see text). 

There, the density of states diverges as in zone II and III. 
The second one arises at Eq = — h where 

J^(-h) = l--^=, (24) 

but, at this energy, the density of states is discontinuous 
in the thermodynamical limit as can be seen in Fig. [2] 

We have confirmed this anomalous behavior numeri- 
cally and observed an even more surprising result. In- 
deed, the energy difference between two consecutive lev- 
els AW — — E^(i = 1, .„, s) is normally given 
by Ao(eo) = l/po(so) in the thermodynamic limit. How- 
ever, in region IV(e), numerical results (at fnite s) show 
that it is not the case. Instead, A" spreads over two 
branches (+) and (— ), depending on the parity of the i. 
In addition, these branches oscillate without converging 
when s increases as can be seen in FigO In this case, the 
gap we computed , in the thermodynamical limit, is the 
average gap, namely A (e ) = ^[^ + \e Q ) + A (_) (e )]. 

To understand physically this unusual phenomenon, 
we have analyzed the classical trajectories in this region 
and have observed that, there are two possible classical 
trajectories [22| ■ Using the results described in Ref. [2(| , 
we computed analytically the expectation values of sev- 
eral observables (such as the magnetization) as a function 
of energy. We found that these values also depend on the 
parity of the level considered but, contrary to the gap, 
the two branches, converge in the thermodynamical limit. 

Finally, note that for h — 0, the LMG model in this re- 
gion coincides with the quantum asymmetric rotor model 
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